# library(devtools)
# install_github('MacoskoLab/liger')
library(tidyverse)
── Attaching packages ───────────────────────────────────────────────────────────────────────────────────── tidyverse 1.2.1 ──
✔ ggplot2 3.2.1     ✔ purrr   0.3.2
✔ tibble  2.1.3     ✔ dplyr   0.8.3
✔ tidyr   1.0.0     ✔ stringr 1.4.0
✔ readr   1.3.1     ✔ forcats 0.4.0
── Conflicts ──────────────────────────────────────────────────────────────────────────────────────── tidyverse_conflicts() ──
✖ tidyr::expand() masks Matrix::expand()
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag()    masks stats::lag()
✖ tidyr::pack()   masks Matrix::pack()
✖ tidyr::unpack() masks Matrix::unpack()
library(liger)
library(Rtsne)
library(reshape2)

Attaching package: ‘reshape2’

The following object is masked from ‘package:tidyr’:

    smiths
library(ggrepel)

Testing LIGER for integration of scATAC and scRNA data. Following vignette from (gitHub)[https://macoskolab.github.io/liger/walkthrough_rna_atac.html].

1. Download and load datasets:

ATAC data is preprocessed to gene-level features (no. of reads that overlap gene body and promoter for each gene). Data is downloaded from (here)[https://umich.app.box.com/s/5hkhpou4inrulo570yhc38ay8g3ont5b]

rna_clusts = readRDS("~/10X_data/rna_cluster_assignments.RDS")
atac_clusts = readRDS("~/10X_data/atac_cluster_assignments.RDS")
pbmc.atac <- readRDS('~/10X_data/pbmc.atac.expression.mat.RDS')
pbmc.rna <- readRDS('~/10X_data/pbmc.rna.expression.mat.RDS')

2. Create LIGER object

List of count matrices as input. createLiger converts to sparseMatrices.

pbmc.data = list(atac=pbmc.atac[,names(atac_clusts)], rna=pbmc.rna[,names(rna_clusts)])
int.pbmc <- createLiger(pbmc.data)
[1] "Removing 97 genes not expressing in atac."
[1] "Removing 11048 genes not expressing in rna."

3. Data preprocessing

## Normalize to column total (adj. for sequencing depth)
int.pbmc <- normalize(int.pbmc)
## Select highly variable genes ONLY IN THE RNA DATASET
int.pbmc <- selectGenes(int.pbmc, datasets.use = 2)
## Scale to the SD
int.pbmc <- scaleNotCenter(int.pbmc)

4. Run NMF

int.pbmc <- optimizeALS(int.pbmc, k=20)

  |                                                                                                                           
  |                                                                                                                     |   0%
  |                                                                                                                           
  |=                                                                                                                    |   1%
  |                                                                                                                           
  |==                                                                                                                   |   2%
  |                                                                                                                           
  |====                                                                                                                 |   3%
  |                                                                                                                           
  |=====                                                                                                                |   4%
  |                                                                                                                           
  |======                                                                                                               |   5%
  |                                                                                                                           
  |=======                                                                                                              |   6%
  |                                                                                                                           
  |========                                                                                                             |   7%
  |                                                                                                                           
  |=========                                                                                                            |   8%
  |                                                                                                                           
  |===========                                                                                                          |   9%
  |                                                                                                                           
  |============                                                                                                         |  10%
  |                                                                                                                           
  |=============                                                                                                        |  11%
  |                                                                                                                           
  |==============                                                                                                       |  12%
  |                                                                                                                           
  |===============                                                                                                      |  13%
  |                                                                                                                           
  |================                                                                                                     |  14%
  |                                                                                                                           
  |==================                                                                                                   |  15%
  |                                                                                                                           
  |===================                                                                                                  |  16%
  |                                                                                                                           
  |=====================================================================================================================| 100%
Converged in 3.941199 mins, 16 iterations.
Best results with seed 1.
smp <- sample(1:nrow(int.pbmc@H$atac), size = 1000)
smp_genes <- sample(1:length(int.pbmc@var.genes), size = 100)
int.pbmc@W[,smp_genes] %>% heatmap(Rowv = NA)
Error in int.pbmc@W[, smp_genes] %>% heatmap(Rowv = NA) : 
  could not find function "%>%"

Align factor

Dodgy quantile normalization

int.pbmc <- readRDS(file = "~/10X_data/pbmc.trainedLIGER.RDS")
cannot open compressed file '/home/jovyan/10X_data/pbmc.trainedLIGER.RDS', probable reason 'No such file or directory'Error in gzfile(file, "rb") : cannot open the connection

int.pbmc@scale.data$atac[,
                         smp_genes] %>%
  melt(varnames=c("cell", "gene")) %>%
  ggplot(aes(value, color=gene)) + geom_histogram() +
  facet_wrap(gene~.)
int.pbmc@scale.data$rna[,
                         smp_genes] %>%
  melt(varnames=c("cell", "gene")) %>%
  ggplot(aes(value, color=gene)) + geom_histogram() +
  facet_wrap(gene~.)

  

Visualize dataset-specific VS common latent factors

Even in dataset specific component there is separation of cell types (Meaningful information that is ATAC specific?)

atac.comp <- int.pbmc@H$atac[smp,] %*% int.pbmc@V$atac
Error: object 'smp' not found
rna.comp <- int.pbmc@H$rna[smp,] %*% int.pbmc@V$rna
rna.common.comp <- int.pbmc@H$rna[smp,] %*% int.pbmc@W
tsne.rna.comp <- Rtsne(rna.comp, pca=F, check_duplicates = F, theta=0.5, perplexity=30)
rownames(tsne.rna.comp$Y) <- rownames(rna.comp)
tsne.rna.common.comp <- Rtsne(rna.common.comp, pca=F, check_duplicates = F, theta=0.5, perplexity=30)
rownames(tsne.rna.common.comp$Y) <- rownames(rna.common.comp)
tsne.rna.comp$Y %>%
  as.tibble(rownames="cell") %>%
  mutate(clust=rna_clusts[cell]) %>%
  ggplot(aes(V1, V2, color=clust)) + 
  geom_point()

tsne.rna.common.comp$Y %>%
  as.tibble(rownames="cell") %>%
  mutate(clust=rna_clusts[cell]) %>%
  ggplot(aes(V1, V2, color=clust)) + 
  geom_point()

Compare weights of dataset specific and common factors

Are HVGs in RNA low coverage in ATAC?

Evaluate feature selection for ATAC-seq data

LS0tCnRpdGxlOiAiVGVzdGluZyBMSUdFUiIKb3V0cHV0OiBodG1sX25vdGVib29rCi0tLQoKYGBge3J9CiMgbGlicmFyeShkZXZ0b29scykKIyBpbnN0YWxsX2dpdGh1YignTWFjb3Nrb0xhYi9saWdlcicpCmxpYnJhcnkodGlkeXZlcnNlKQpsaWJyYXJ5KGxpZ2VyKQpsaWJyYXJ5KFJ0c25lKQpsaWJyYXJ5KHJlc2hhcGUyKQpsaWJyYXJ5KGdncmVwZWwpCmBgYAoKVGVzdGluZyBMSUdFUiBmb3IgaW50ZWdyYXRpb24gb2Ygc2NBVEFDIGFuZCBzY1JOQSBkYXRhLiBGb2xsb3dpbmcgdmlnbmV0dGUgZnJvbSAoZ2l0SHViKVtodHRwczovL21hY29za29sYWIuZ2l0aHViLmlvL2xpZ2VyL3dhbGt0aHJvdWdoX3JuYV9hdGFjLmh0bWxdLgoKIyMjIDEuIERvd25sb2FkIGFuZCBsb2FkIGRhdGFzZXRzOiAKQVRBQyBkYXRhIGlzIHByZXByb2Nlc3NlZCB0byBnZW5lLWxldmVsIGZlYXR1cmVzIChuby4gb2YgcmVhZHMgdGhhdCBvdmVybGFwIGdlbmUgYm9keSBhbmQgcHJvbW90ZXIgZm9yIGVhY2ggZ2VuZSkuIERhdGEgaXMgZG93bmxvYWRlZCBmcm9tIChoZXJlKVtodHRwczovL3VtaWNoLmFwcC5ib3guY29tL3MvNWhraHBvdTRpbnJ1bG81NzB5aGMzOGF5OGczb250NWJdCmBgYHtyfQpybmFfY2x1c3RzID0gcmVhZFJEUygifi8xMFhfZGF0YS9ybmFfY2x1c3Rlcl9hc3NpZ25tZW50cy5SRFMiKQphdGFjX2NsdXN0cyA9IHJlYWRSRFMoIn4vMTBYX2RhdGEvYXRhY19jbHVzdGVyX2Fzc2lnbm1lbnRzLlJEUyIpCnBibWMuYXRhYyA8LSByZWFkUkRTKCd+LzEwWF9kYXRhL3BibWMuYXRhYy5leHByZXNzaW9uLm1hdC5SRFMnKQpwYm1jLnJuYSA8LSByZWFkUkRTKCd+LzEwWF9kYXRhL3BibWMucm5hLmV4cHJlc3Npb24ubWF0LlJEUycpCmBgYAoKIyMjIDIuIENyZWF0ZSBMSUdFUiBvYmplY3QKTGlzdCBvZiBjb3VudCBtYXRyaWNlcyBhcyBpbnB1dC4gYGNyZWF0ZUxpZ2VyYCBjb252ZXJ0cyB0byBzcGFyc2VNYXRyaWNlcy4gCmBgYHtyfQpwYm1jLmRhdGEgPSBsaXN0KGF0YWM9cGJtYy5hdGFjWyxuYW1lcyhhdGFjX2NsdXN0cyldLCBybmE9cGJtYy5ybmFbLG5hbWVzKHJuYV9jbHVzdHMpXSkKaW50LnBibWMgPC0gY3JlYXRlTGlnZXIocGJtYy5kYXRhKQoKYGBgCiMjIyAzLiBEYXRhIHByZXByb2Nlc3NpbmcKYGBge3J9CiMjIE5vcm1hbGl6ZSB0byBjb2x1bW4gdG90YWwgKGFkai4gZm9yIHNlcXVlbmNpbmcgZGVwdGgpCmludC5wYm1jIDwtIG5vcm1hbGl6ZShpbnQucGJtYykKIyMgU2VsZWN0IGhpZ2hseSB2YXJpYWJsZSBnZW5lcyBPTkxZIElOIFRIRSBSTkEgREFUQVNFVAppbnQucGJtYyA8LSBzZWxlY3RHZW5lcyhpbnQucGJtYywgZGF0YXNldHMudXNlID0gMikKIyMgU2NhbGUgdG8gdGhlIFNECmludC5wYm1jIDwtIHNjYWxlTm90Q2VudGVyKGludC5wYm1jKQpgYGAKCiMjIyA0LiBSdW4gTk1GIAoKYGBge3J9CmludC5wYm1jIDwtIG9wdGltaXplQUxTKGludC5wYm1jLCBrPTIwKQpgYGAKYGBge3J9CnNtcCA8LSBzYW1wbGUoMTpucm93KGludC5wYm1jQEgkYXRhYyksIHNpemUgPSAxMDAwKQpzbXBfZ2VuZXMgPC0gc2FtcGxlKDE6bGVuZ3RoKGludC5wYm1jQHZhci5nZW5lcyksIHNpemUgPSAxMDApCgppbnQucGJtY0BXWyxzbXBfZ2VuZXNdICU+JSBoZWF0bWFwKFJvd3YgPSBOQSkKaW50LnBibWNAViRhdGFjWyxzbXBfZ2VuZXNdICU+JSBoZWF0bWFwKFJvd3YgPSBOQSkKYGBgCmBgYHtyfQoKYGBgCgojIyMgQWxpZ24gZmFjdG9yCkRvZGd5IHF1YW50aWxlIG5vcm1hbGl6YXRpb24KYGBge3J9CmludC5wYm1jIDwtIHF1YW50aWxlQWxpZ25TTkYoaW50LnBibWMpCiMgc2F2ZVJEUyhpbnQucGJtYywgZmlsZSA9ICJ+LzEwWF9kYXRhL3BibWMudHJhaW5lZExJR0VSLlJEUyIpCmludC5wYm1jIDwtIHJlYWRSRFMoZmlsZSA9ICJ+L215X2RhdGEvMTBYX2RhdGEvcGJtYy50cmFpbmVkTElHRVIuUkRTIikKYGBgCgo8IS0tICMjIyBWaXN1YWxpemUgIC0tPgo8IS0tIGBgYHtyfSAtLT4KPCEtLSBydW5NeVRTTkUgPC0gZnVuY3Rpb24ob2JqZWN0LCB1c2UucmF3ID0gRiwgZGltcy51c2UgPSAxOm5jb2wob2JqZWN0QEgubm9ybSksIHVzZS5wY2EgPSBGLCAtLT4KPCEtLSAgICAgICAgICAgICAgICAgICAgIHBlcnBsZXhpdHkgPSAzMCwgdGhldGEgPSAwLjUsIG1ldGhvZCA9ICdSdHNuZScsIGZpdHNuZS5wYXRoID0gTlVMTCwgLS0+CjwhLS0gICAgICAgICAgICAgICAgICAgICByYW5kLnNlZWQgPSA0MikgeyAtLT4KPCEtLSAgIGRhdGEudXNlIDwtIGRvLmNhbGwocmJpbmQsIG9iamVjdEBIKSAgIC0tPgo8IS0tICAgc2V0LnNlZWQocmFuZC5zZWVkKSAtLT4KPCEtLSAgIG9iamVjdEB0c25lLmNvb3JkcyA8LSBSdHNuZShvYmplY3RASC5ub3JtWywgZGltcy51c2VdLCBwY2EgPSB1c2UucGNhLCBjaGVja19kdXBsaWNhdGVzID0gRiwgLS0+CjwhLS0gICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICB0aGV0YSA9IHRoZXRhLCBwZXJwbGV4aXR5ID0gcGVycGxleGl0eSkkWSAtLT4KPCEtLSAgIHJvd25hbWVzKG9iamVjdEB0c25lLmNvb3JkcykgPC0gcm93bmFtZXMoZGF0YS51c2UpIC0tPgo8IS0tICAgcmV0dXJuKG9iamVjdCkgLS0+CjwhLS0gfSAtLT4KCmBgYHtyfQpzbXBfZ2VuZXMgPC0gc2FtcGxlKDE6MjAwMCwxMCkKaW50LnBibWNAc2NhbGUuZGF0YSRhdGFjWywKICAgICAgICAgICAgICAgICAgICAgICAgIHNtcF9nZW5lc10gJT4lCiAgbWVsdCh2YXJuYW1lcz1jKCJjZWxsIiwgImdlbmUiKSkgJT4lCiAgZ2dwbG90KGFlcyh2YWx1ZSwgY29sb3I9Z2VuZSkpICsgZ2VvbV9oaXN0b2dyYW0oKSArCiAgZmFjZXRfd3JhcChnZW5lfi4pCgppbnQucGJtY0BzY2FsZS5kYXRhJHJuYVssCiAgICAgICAgICAgICAgICAgICAgICAgICBzbXBfZ2VuZXNdICU+JQogIG1lbHQodmFybmFtZXM9YygiY2VsbCIsICJnZW5lIikpICU+JQogIGdncGxvdChhZXModmFsdWUsIGNvbG9yPWdlbmUpKSArIGdlb21faGlzdG9ncmFtKCkgKwogIGZhY2V0X3dyYXAoZ2VuZX4uKQogIApgYGAKCgo8IS0tIGBgYCAtLT4KCiMjIyBWaXN1YWxpemUgZGF0YXNldC1zcGVjaWZpYyBWUyBjb21tb24gbGF0ZW50IGZhY3RvcnMKCkV2ZW4gaW4gZGF0YXNldCBzcGVjaWZpYyBjb21wb25lbnQgdGhlcmUgaXMgc2VwYXJhdGlvbiBvZiBjZWxsIHR5cGVzIChNZWFuaW5nZnVsIGluZm9ybWF0aW9uIHRoYXQgaXMgQVRBQyBzcGVjaWZpYz8pCgpgYGB7cn0KCmF0YWMuY29tcCA8LSBpbnQucGJtY0BIJGF0YWNbc21wLF0gJSolIGludC5wYm1jQFYkYXRhYwphdGFjLmNvbW1vbi5jb21wIDwtIGludC5wYm1jQEgkYXRhY1tzbXAsXSAlKiUgaW50LnBibWNAVwoKdHNuZS5hdGFjLmNvbXAgPC0gUnRzbmUoYXRhYy5jb21wLCBwY2E9VCwgY2hlY2tfZHVwbGljYXRlcyA9IEYsIHRoZXRhPTAuNSwgcGVycGxleGl0eT0xMCkKcm93bmFtZXModHNuZS5hdGFjLmNvbXAkWSkgPC0gcm93bmFtZXMoYXRhYy5jb21wKQp0c25lLmF0YWMuY29tbW9uLmNvbXAgPC0gUnRzbmUoYXRhYy5jb21tb24uY29tcCwgcGNhPVQsIGNoZWNrX2R1cGxpY2F0ZXMgPSBGLCB0aGV0YT0wLjUsIHBlcnBsZXhpdHk9MzApCnJvd25hbWVzKHRzbmUuYXRhYy5jb21tb24uY29tcCRZKSA8LSByb3duYW1lcyhhdGFjLmNvbW1vbi5jb21wKQoKCnRzbmUuYXRhYy5jb21wJFkgJT4lCiAgYXMudGliYmxlKHJvd25hbWVzPSJjZWxsIikgJT4lCiAgbXV0YXRlKGNsdXN0PWF0YWNfY2x1c3RzW2NlbGxdKSAlPiUKICBnZ3Bsb3QoYWVzKFYxLCBWMiwgY29sb3I9Y2x1c3QpKSArIAogIGdlb21fcG9pbnQoKQoKdHNuZS5hdGFjLmNvbW1vbi5jb21wJFkgJT4lCiAgYXMudGliYmxlKHJvd25hbWVzPSJjZWxsIikgJT4lCiAgbXV0YXRlKGNsdXN0PWF0YWNfY2x1c3RzW2NlbGxdKSAlPiUKICBnZ3Bsb3QoYWVzKFYxLCBWMiwgY29sb3I9Y2x1c3QpKSArIAogIGdlb21fcG9pbnQoKQoKYGBgCgpgYGB7cn0KCnJuYS5jb21wIDwtIGludC5wYm1jQEgkcm5hW3NtcCxdICUqJSBpbnQucGJtY0BWJHJuYQpybmEuY29tbW9uLmNvbXAgPC0gaW50LnBibWNASCRybmFbc21wLF0gJSolIGludC5wYm1jQFcKCnRzbmUucm5hLmNvbXAgPC0gUnRzbmUocm5hLmNvbXAsIHBjYT1GLCBjaGVja19kdXBsaWNhdGVzID0gRiwgdGhldGE9MC41LCBwZXJwbGV4aXR5PTMwKQpyb3duYW1lcyh0c25lLnJuYS5jb21wJFkpIDwtIHJvd25hbWVzKHJuYS5jb21wKQp0c25lLnJuYS5jb21tb24uY29tcCA8LSBSdHNuZShybmEuY29tbW9uLmNvbXAsIHBjYT1GLCBjaGVja19kdXBsaWNhdGVzID0gRiwgdGhldGE9MC41LCBwZXJwbGV4aXR5PTMwKQpyb3duYW1lcyh0c25lLnJuYS5jb21tb24uY29tcCRZKSA8LSByb3duYW1lcyhybmEuY29tbW9uLmNvbXApCgoKdHNuZS5ybmEuY29tcCRZICU+JQogIGFzLnRpYmJsZShyb3duYW1lcz0iY2VsbCIpICU+JQogIG11dGF0ZShjbHVzdD1ybmFfY2x1c3RzW2NlbGxdKSAlPiUKICBnZ3Bsb3QoYWVzKFYxLCBWMiwgY29sb3I9Y2x1c3QpKSArIAogIGdlb21fcG9pbnQoKQoKdHNuZS5ybmEuY29tbW9uLmNvbXAkWSAlPiUKICBhcy50aWJibGUocm93bmFtZXM9ImNlbGwiKSAlPiUKICBtdXRhdGUoY2x1c3Q9cm5hX2NsdXN0c1tjZWxsXSkgJT4lCiAgZ2dwbG90KGFlcyhWMSwgVjIsIGNvbG9yPWNsdXN0KSkgKyAKICBnZW9tX3BvaW50KCkKCmBgYAoKIyMjIENvbXBhcmUgd2VpZ2h0cyBvZiBkYXRhc2V0IHNwZWNpZmljIGFuZCBjb21tb24gZmFjdG9ycwpgYGB7cn0KaW50LnBibWNAVyAlPiUgCiAgbWVsdCh2YXJuYW1lcyA9IGMoImZhY3RvciIsICJnZW5lIikpICU+JQogIGZpbHRlcihmYWN0b3I9PTEpICU+JQogIG11dGF0ZShyYW5rPXJhbmsodmFsdWUpKSAlPiUKICBnZ3Bsb3QoYWVzKHJhbmssIHZhbHVlKSkgKyAKICBnZW9tX3BvaW50KCkgKwogIGdlb21fdGV4dF9yZXBlbChkYXRhPS4gJT4lIGZpbHRlcihyYW5rID4gMjMwMCksIGFlcyhsYWJlbD1nZW5lKSkKYGBgCmBgYHtyfQpWX2F0YWMgPC0gaW50LnBibWNAViRhdGFjIApjb2xuYW1lcyhWX2F0YWMpIDwtIGNvbG5hbWVzKGludC5wYm1jQFcpCgpWX3JuYSA8LSBpbnQucGJtY0BWJHJuYQpjb2xuYW1lcyhWX3JuYSkgPC0gY29sbmFtZXMoaW50LnBibWNAVykKClZfYXRhYyAlPiUgCiAgbWVsdCh2YXJuYW1lcyA9IGMoImZhY3RvciIsICJnZW5lIikpICU+JQogIGZpbHRlcihmYWN0b3I9PTEpICU+JQogIG11dGF0ZShyYW5rPXJhbmsodmFsdWUpKSAlPiUKICBnZ3Bsb3QoYWVzKHJhbmssIHZhbHVlKSkgKyAKICBnZW9tX3BvaW50KCkgKwogIGdlb21fdGV4dF9yZXBlbChkYXRhPS4gJT4lIGZpbHRlcihyYW5rID4gMjMwMCksIGFlcyhsYWJlbD1nZW5lKSkgKwogIGdndGl0bGUoJ3RvcCBWIEFUQUMnKQoKVl9ybmEgJT4lIAogIG1lbHQodmFybmFtZXMgPSBjKCJmYWN0b3IiLCAiZ2VuZSIpKSAlPiUKICBncm91cF9ieShmYWN0b3IpICU+JQogIG11dGF0ZShyYW5rPXJhbmsodmFsdWUpKSAlPiUKICB1bmdyb3VwKCkgJT4lCiAgZ2dwbG90KGFlcyhyYW5rLCB2YWx1ZSkpICsgCiAgZ2VvbV9wb2ludCgpICsKICAjIGdlb21fdGV4dF9yZXBlbChkYXRhPS4gJT4lIGZpbHRlcihyYW5rID4gMjMxMCksIGFlcyhsYWJlbD1nZW5lKSkgKwogIGZhY2V0X3dyYXAoZmFjdG9yfi4pICsKICBnZ3RpdGxlKCd0b3AgViBBVEFDJykKCmBgYApgYGB7cn0KZ2VuZXMuYmMgPC0gcmVhZC50YWJsZShmaWxlID0gIi4uL215X2RhdGEvY2VsbHJhbmdlci1hdGFjMTEwX2NvdW50XzMwNDM5X1dTU1M4MDM4MzYwX0dSQ2gzOC0xXzFfMC5nZW5lc19iYy5iZWQiLCBzZXAgPSAiXHQiLCBhcy5pcyA9IGMoMSksIGhlYWRlciA9IEZBTFNFKQoKZ2VuZXMuYmMKYGBgCgojIyMgQXJlIEhWR3MgaW4gUk5BIGxvdyBjb3ZlcmFnZSBpbiBBVEFDPwpgYGB7cn0KCmxvbmdfYXRhY19yYXcgPC0gCiAgaW50LnBibWNAcmF3LmRhdGEkYXRhYyAlPiUKICBhcy5tYXRyaXgoKSAlPiUKICBtZWx0KHZhcm5hbWVzPWMoImdlbmUiLCAiY2VsbCIpKSAlPiUKICBtdXRhdGUodmFyX2dlbmU9aWZlbHNlKGdlbmUgJWluJSBpbnQucGJtY0B2YXIuZ2VuZXMsIFQsIEYpKSAKCmNlbGxfc21wIDwtIHNhbXBsZSh1bmlxdWUobG9uZ19hdGFjX3JhdyRjZWxsKSwgNTAwKQpnZW5lX3NtcCA8LSBzYW1wbGUodW5pcXVlKGxvbmdfYXRhY19yYXckZ2VuZSksIDEwMDApCgpsb25nX2F0YWNfcmF3ICU+JQogIGZpbHRlcihnZW5lICVpbiUgZ2VuZV9zbXAgJiBjZWxsICVpbiUgY2VsbF9zbXApICU+JQogIGdncGxvdChhZXModmFsdWUsIGZpbGw9dmFyX2dlbmUpKSArIAogIGdlb21faGlzdG9ncmFtKGJpbnM9MTAwKSArCiAgZmFjZXRfd3JhcCh2YXJfZ2VuZX4uLCBzY2FsZXMgPSAiZnJlZV95IikgKwogIGNvb3JkX2NhcnRlc2lhbih4bGltPWMoMCwzMCkpCgoKYGBgCgpgYGB7cn0KbG9uZ19hdGFjX3JhdyAlPiUKICBmaWx0ZXIoZ2VuZSAlaW4lIGdlbmVfc21wICYgY2VsbCAlaW4lIGNlbGxfc21wKSAlPiUKICBtdXRhdGUoemVyb2VzPWlmZWxzZSh2YWx1ZT09MCwgVCxGKSkgJT4lCiAgZ2dwbG90KGFlcyh2YXJfZ2VuZSwgZmlsbD16ZXJvZXMpKSArIAogIGdlb21fYmFyKHBvc2l0aW9uPSJmaWxsIikKYGBgCgoKIyMjIEV2YWx1YXRlIGZlYXR1cmUgc2VsZWN0aW9uIGZvciBBVEFDLXNlcSBkYXRhCmBgYHtyfQppbnQucGJtY19odmdfcm5hIDwtIHNlbGVjdEdlbmVzKGludC5wYm1jLCBudW0uZ2VuZXMgPSAyMDAwLCBkYXRhc2V0cy51c2UgPSAyKQppbnQucGJtY19odmdfYXRhYyA8LSBzZWxlY3RHZW5lcyhpbnQucGJtYywgbnVtLmdlbmVzID0gMjAwMCwgZGF0YXNldHMudXNlID0gMSkKCnBsb3RfbWVhbkFjY19odmdzIDwtIGZ1bmN0aW9uKGludC5wYm1jKXsKICBnZW5lX21lYW5zIDwtIAogICAgaW50LnBibWNAcmF3LmRhdGEkYXRhYyAlPiUKICAgIHJvd01lYW5zKCkKICBkYXRhLmZyYW1lKG1lYW5BY2M9Z2VuZV9tZWFucywgZ2VuZT1uYW1lcyhnZW5lX21lYW5zKSkgJT4lCiAgICAgbXV0YXRlKHZhcl9nZW5lPWlmZWxzZShnZW5lICVpbiUgaW50LnBibWNAdmFyLmdlbmVzLCAiSFZHIiwgIm90aGVyIikpICU+JQogICAgZ2dwbG90KGFlcyhmaWxsPXZhcl9nZW5lLCBtZWFuQWNjKSkgKwogICAgZ2VvbV9oaXN0b2dyYW0oYmlucz01MCkgKwogICAgZmFjZXRfd3JhcCh2YXJfZ2VuZX4uLCBzY2FsZXM9ImZyZWVfeSIsIG5yb3c9MiwgbmNvbD0xKQogIH0KCmdnYXJyYW5nZSgKICBwbG90X21lYW5BY2NfaHZncyhpbnQucGJtY19odmdfcm5hKSArIGdndGl0bGUoIkhWRyBSTkEiKSwKICBwbG90X21lYW5BY2NfaHZncyhpbnQucGJtY19odmdfYXRhYykgKyBnZ3RpdGxlKCJIVkcgQVRBQyIpLCBjb21tb24ubGVnZW5kID0gVFJVRQogICkKYGBgCgoKCgoKCgoK